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We report on theoretical study of second harmonic generation in graphene. Phenomenological 
analysis based on symmetry arguments is carried out. It is demontrated, that in ideal graphene 
samples second harmonic generation is possible only if the radiation wave vector or its magnetic field 
is taken into account. Microscopic theory is developed for the classical regime of radiation interaction 
with electrons, where photon energy is much smaller than the charge carriers characteristic energy. 
It is demonstrated, that the emitted radiation can be strongly circularly polarized for the linearly 
polarized incident wave. 
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Introduction. Graphene, a monolayer of carbon atoms 
arranged in honeycomb lattice, became recently one of 
the most promising materials owing to its intriguing elec- 
tronic and transport properties, for recent reviews see 
Refs. pH!] and references wherein. 

It is well known that non-linear transport and optical 
phenomena are most sensitive to the fine details of en- 
ergy spectrum and carrier scattering in condensed mat- 
ter [5J [5] • One of the non- linear effects is the generation 
of dc current under homogeneous excitation caused by 
the asymmetry of optical transitions and scattering in 
non-centrosymmetric media (photogalvanic effect) or by 
the transfer of photon momentum to the electron sys- 
tem (photon drag effect) [7]. The latter one does not 
require low symmetry of the structure, but it is masked 
by various photogalvanic effects in ordinary semiconduc- 
tor quantum wells which lack an inversion center as a 
rule [5]. By contrast, graphene is a centrosymmetric ma- 
terial which allowed unambiguous identification of the 
photon drag effect by its remarkable polarization depen- 
dence [5]. 

Here we consider theoretically light wave vector in- 
duced second harmonic generation in graphene which 
consists in the formation of the ac response: electric cur- 
rent or dielectric polarization oscillating at the double 
frequency as compared with the frequency of the incident 
radiation. Experimentally second harmonic in graphene 
was observed in [101 lllj , but microscopic theory is absent 
so far. Alike the photon drag effect, in centrosymmetric 
media like graphene second harmonic generation is al- 
lowed only provided that the photon momentum trans- 
fer to electron system is taken into account. This is in 
contrast to the third (and other odd) harmonic genera- 
tion [T21 [13] as well as coherent injection of ballistic cur- 
rents |14[ 1 1 5 j which do not require photon momentum 
transfer, but appear in higher orders in incident electro- 
magnetic field. 

The second harmonic response is calculated here for 
low frequencies, e.g. terahertz radiation. In such case, 
the photon energy is much smaller than characteristic en- 
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ergy of electrons and a classical regime is realized for the 
radiation interaction with electrons in graphene, which 
brings about very strong polarization effects. In particu- 
lar, second harmonic can be strongly circularly polarized 
for linearly polarized incident radiation. 

Phenomenological theory. We consider an ideal single 
layer graphene described by the point symmetry group 
D 6 h. This point group contains an inversion center, 
therefore, the second harmonic generation is possible only 
if one takes into consideration the light wave vector q 
transfer to electron ensemble or joint action of electric 
E and magnetic fields B of the incident radiation [SJ [7J . 
Since the magnetic field components B a are proportional 
to certain combinations qpE 1 , where greek subscripts 
enumerate Cartesian coordinates, in phenomenological 
considerations it is enough to allow for the contributions 
to the second order response being linear in the wave vec- 
tor components and quadratic in the electric field ampli- 
tudes. 

It is convenient to study the electric current, j(2u>), in- 
duced at the double frequency 2uj as a response to the in- 
cident radiation with the fundamental frequency uj. The 
electric current is related with the oscillating polariza- 
tion, P(2cj), as j(2oj) = —2\u)P{2ui). The general phe- 
nomenological expression describing the second harmonic 
generation reads 

j\(2u) = SxaunQuEuE^. (1) 

Here Sx^rj is the material tensor and the electric field of 
the fundamental wave is taken in the complex form 

E(r,t) = Ee Kqr ~ ut) + c.c. . (2) 

The same complex form, Eq. ([2]), is assumed for all other 
quantities in question: the magnetic field, B, and the 
current density, j. It is worth mentioning that, in con- 
trast to the photon drag effect where the dc current is 
proportional to the bilinear combinations of electric field 
and its conjugate, E„E*, the second harmonic generation 
is described by the quadratic combinations E V E V . Hence, 
tensor S\^ vv is invariant with respect to permutations of 
the two last subscripts. The elements S\^ vr , are, in gen- 
eral, complex, the presence of both real and imaginary 
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FIG. 1: Schematic illustration of the second harmonic gener- 
ation geometry, Eqs. Q. It is seen that there are both longi- 
tudinal [panel (a)] and transverse [panel (b)] components of 
the second order response. 



parts describes the phase shift between the fundamental 
and second harmonics. 

There are 7 independent parameters, Si ... S^, which 
describe second harmonic generation in the system of 
point group D^: 

j x {2u) = S x q x [El + E 2 y ) + S 2 [q x (E 2 x - E 2 y ) + 2q y E x E y ] 



S 3 q z E x E z + S 4 q x E 2 z , (3a) 



j y (2u) = S igy (E 2 x + E 2 y ) + S 2 [q y (E 2 - E 2 ) + 2q x E x E y ] 
+ S 3 q z E y E z +S iqy E 2 , (3b) 



j,(2w) = S 5 q z E 2 +S 6 q z (E 2 x 



-S7(q x E x E z +q y E y E z ), 
( 3c ) 

where we assumed that the graphene layer occupies (xy) 
plane. Expressions ( 3a I and ( 3b I are similar to phe- 



nomenological equations describing linear drag effect in 
graphene [16] . 

The specifics of the graphene band structure imposes 
strong restrictions on values of the constants S± , . . . , S7 
in Eqs. (J3|. Indeed, conduction and valence band of 
graphene are formed from 7r orbitals of carbon atoms. 
As a result, the transitions in z polarization are possible 
only with allowance for other, <j, orbitals which are sepa- 
rated from 7T orbitals by several electronvolts. It implies 
that in the model where only tt orbitals are taken into 
account (e.g. in the classical frequency range of interest) 
constants Si = for i > 2 and second harmonic genera- 
tion is described by the simplified expressions with only 
two independent parameters: 



j x (2u) = Siq x (E 2 + E 2 ) 



S 2 [q x (E 2 -E 2 )+2q y E x E y ], 
(4a) 



j v (2u) = Siq y (E 2 + E 2 ) + S 2 [q y (E 2 y - E 2 ) + 2q x E x E y }. 

(4b) 

Figure [l] schematically shows the geometry of the sec- 
ond harmonic generation and the response at a dou- 
ble frequency, 2oj. For instance, for q x 7^ 0, q y = 0, 
there is a component of the in-plane oscillating current 
j(2u>) parallel to the light incidence plane described by 
(S 1 +S 2 )q x E 2 x + (Si-S 2 )q x E 2 [Figjlfa)]. Additionally, 
there is a transverse contribution, 2S 2 q x E x E y , Fig. Qb). 
Microscopic theory for these constants Si and £2 is de- 
veloped below. 

Graphene samples are usually deposited on substrates, 
hence z — > —z symmetry is broken. As a result, structure 



becomes non-centrosymmetric, described by the point 
group Cq v , and second harmonic generation becomes pos- 
sible for q = 0. It can be shown that additional contri- 
butions have the following form 



j x (2u) = S' 3 E X E Z , j y (2u>) 



j z (2u) = S' 5 E 2 z + S' e (E 2 x 



S' 3 E y E z 



(5a) 



(5b) 



where S' 3 , S' 5 and S' 6 are some coefficients directly pro- 
portional to parameters d{ 5 , d 31 and 9g 3 introduced in 
Eq. (10) of Ref. [TT]. As one can see from Eqs. ^ these 
contributions require transitions in z polarization and are 
strongly suppressed in classical frequency range [16] but 
they may take over for higher frequencies [llj . 

Microscopic theory. We consider here a doped 
graphene layer and assume that the characteristic elec- 
tron energy (being the Fermi energy, sf, for degenerate 
and the temperature, k^T, for non-degenerate electrons) 
is much larger than the photon energy hu. This condi- 
tion allows one to treat electron dynamics in the frame- 
work of Boltzmann equation. For typical carrier densi- 
ties on the order of 10 12 cm -2 the electron Fermi energy 
£p ~ 100 meV, therefore carriers are degenerate even 
at room temperature and classical treatment is valid for 
microwave or terahertz radiation [9]. 

We introduce time, t, position, r, and momentum, p, 
dependent electron distribution function f(p,r, t) which 
satisfies standard kinetic equation US] 



df 

at 



' dr 



-[v p x B] 

c 



df 
dp 



Q{f}- (6) 



Here E and B are the electric and magnetic fields of the 
incident wave taken in the form Eq. ([2]), e is the electron 
charge, e = — |e|, c is the speed of light, v p is the electron 
velocity in the given state p, 



dp 



v P 



\P\ 



(7) 



Here e p = v\p\ is the electron dispersion, v « c/300 is 
effective electron speed in graphene, and Q{f} in Eq. Q 
is the collision integral. In Eq. ^ it is assumed that 
f(p,r, t) describes the electron distribution in any of two 
equivalent degenerate valleys, K and K'. It is clearly 
seen from Eq. that electron motion is affected by the 
in-plane components of electric field E and by z compo- 
nent of magnetic field B. 

The distribution function is sought in a form |16] 

f(k,r,t) = Me) + [f^(k)e^ r -^ +c.c.] + 

/ 2 °(fe) + [/ 2 2 -(fc) e 2i ^-^+c.c] , (8) 

where /o(e) is the equilibrium (Fcrmi-Dirac) distribution 
function, /"(fe) is the first order correction, and functions 
/ 2 °, /|" are the second order corrections corresponding to 
the responses at zero frequency (photon drag effect) and 
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at the frequency 2ui (second harmonic) , respectively. The 
terms containing higher frequencies (and, hence, higher 
powers of E, B) are disregarded here because they do 
not contribute to the second harmonic generation. The 
dc response was calculated in Refs. [TS] , below we find 
the function /f w (fc). 

In order to determine /| w (fe) we solve Eq. (|6| by iter- 
ations. To that end one has to find first order correction, 
fx(k), which can be expressed up to the linear in q terms 

as nsi 



fx(k) = -er ltU fo\ (E v p ) 



[(qv p )(E v p ) - v (qE )/2] 



v 2 {qE ) 
2uj 



(9) 



Here we introduced f = d/o/de, 



1 



(n = l,2), 



where T\ and r 2 are the relaxation times for the first and 
second angular harmonics of the distribution function 
which describe the momentum (first angular harmonic) 
and momentum alignment (second angular harmonic) re- 
laxation of electrons, respectively. These relaxation times 
can be determined in Born approximation as |17j 



1 



1 

T n {e k ) 



2tt x - 

-^-"■imp 

P 



\Vh-t 



1 + COS 7? 



(10) 

[1 - cos (ni9)]6(e k - e p ), 



where rii mp is the (two-dimensional) impurity density, V q 
is the Fourier-transform of the impurity potential, ■& is 
the scattering angle, (l + cos$)/2 comes from the overlap 
of Bloch amplitudes and the normalization area is set 
to unity hereafter. Intervalley scattering processes are 
disregarded in Eq. (10 1. 



The energy relaxation processes are neglected assum- 
ing that ujT e ^ 1, Te/ti^ 3> 1, where t £ is the en- 
ergy relaxation time. Equation ^ is valid as long as 
qvTi, qvT2, qv juj <C 1, that is if spatial variations of the 
distribution function occur on the scale exceeding by far 
the mean free path, this condition is readily fulfilled in 
the relevant frequency range. 

The second order correction satisfies the following 
equation 

-2iujfi"(p) + 2iqv p fl'»(p)+ 

et+l[ flp xB])^=Q{n (11) 

and the electric current density is given by: 

j=4e£vf"(p), (12) 



(a) 






0.0 


■a 


-0.2 




-0.4 


I 


-0.6 




-0.8 








-1 n 









1 / ' ,- 

\ 1 








\ 1 / 







































(b) l.Op 

i 0.6 

4 0.4 

5 0.2 1 
6? 0.0 
* -0.2- 



1 2 3 4 5 
Radiation Frequency (wri) 




















































































2 3 4 5 



Radiation Frequency (ojti) 



FIG. 2: Frequency dependence of the contributions to the 
second harmonic response: qSi is shown in panel (a), and qS2 
is shown in panel (b). Solid blue lines depict real part and 
dashed red lines depict imaginary part of the response. Values 
of current are given in arbitrary units, radiation frequency is 
presented in units of ojti, where ri is momentum relaxation 
time. Calculation is carried out for the short-range scattering 
where n = 2t2 oc e^T 1 for degenerate electrons. 



where the factor 4 takes into account spin and valley 
degeneracy. Solving Eq. (Ill, retaining in f 2 u {p) only 



linear in q or linear in B terms and expressing B = 
[(q/l) x E] (f° r the sample in vacuum) we arrive at the 
following expressions for the constants S\ and 62: 



Si 



e 3 v 4 
2oj 



Jo* 



(13a) 



T\2u 



S 2 = 



2cu 



(3 + ilut 2 m) + (1 - iwT 2lUJ ) 



dc-fe 



(13b) 



-^-(1 + 4\lot 2 ,2lj) - ^-[ri,2 W (l - 2iwr 2 .2 W )] 
e k de k 

It is noteworthy, that the physical mechanisms govern- 
ing second harmonic generation are the same as in the 
photon drag effect [TH]. First one is related to the joint 
action of electric and magnetic fields of the incident ra- 
diation: electric field induces ac current with the funda- 
mental frequency u), while magnetic field oscillating at 
the same frequency induces both dc current and the cur- 
rent component with the double frequency, 2u>. This is 
EB mechanism of the second harmonic generation. An- 
other qE 2 mechanism is connected with the spatial de- 
pendence of the electric field, its qualitative description 
can be carried out along the same lines as in Ref. |16j . 
In quantum language the microscopic mechanisms gov- 
erning these effects are due to the magneto-dipole and 
quadrupole transitions, respectively. 

Figure [2] shows frequency dependence of electric cur- 
rent components at double frequency j(2u>) proportional 
to coefficients Si [panel (a)] and S2 [panel (b)], respec- 
tively. Calculation is carried out for the short-range scat- 
tering where Vk- P in Eq. ( 10 1 is constant and relaxation 



times arc interrelated as t\ = 2r 2 cx . The elec- 
tron gas is assumed to be degenerate. It can be seen 
that the response at double frequency has complex non- 
monotonous behavior. In the static limit, uj —¥ 0, the co- 
efficients Si and S2 are real and diverge as 1/w, but the 
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net current remains finite due to factors oc q in Eqs. @. 
Coefficients Si and S2 become, up to common factor, 
equal to the constants Ti and T 2 describing linear photon 
drag effect, see Eqs. (26) of Ref. (T5J, because at lo = re- 
sponses at zero and double frequencies arc indistinguish- 
able. At high frequencies, ujt± 3> 1, wt 2 3> 1, parameters 
Si and are proportional to 1/w 3 , hence, current den- 
sity decays as 1/ui 2 . The imaginary part of current is 
determined by the retardation between the driving elec- 
tric held and electron velocity [9j[T6]. Therefore, it is zero 
in the limit of static helds, shows maximum at u)T\ ~ 1 
and decreases for cjti » 1. Typical values of j(2u>) for 
luti ~ 1 are close to those of dc current. For the experi- 
mental parameters of Ref. [9]: j(2u))/P ~ 0.5 nAcm/W, 
where P is the power density. 

Polarization of the second harmonic. The electromag- 
netic field emitted at the double frequency can be found 
by means of Maxwell equations. It is convenient to deter- 
mine the vector potential A(r,t) oscillating at 2cj, which 
satisfies the following equation [T8] 



4w 2 

AA{r,t) + — A{2u) 



4tt 



,2it/|i p— 2iujt 



6(z)j<M 



(14) 

where q\\ — (q x ,q y ) is the projection of the incident radia- 
tion wave vector onto the graphene layer and we assumed 
that the sample is positioned in vacuum at z = plane. 
Using one-dimensional Green's function [7], the solution 
of Eq. ( 14 1 can be written as 



A(r,t) = — <> 



2iqf|| p—litot 

cq. " 
3(2")q\i , >„, 



'! 



j(2w)e 2i ^ |z| 
(e^ |z| - l)signz + c.c. 



Here e z is the unit vector along z-axis. It follows from 
Eq. (15) that the emitted radiation propagates with the 



in-plane wave vector 2qy . The dependence of the vec- 
tor potential on z given by the factor exp (2iq^|z|) corre- 
sponds to the wave diverging from the quantum well. 

The electric field of the emitted second harmonic is re- 
lated with its vector potential by the standard expression 
E(2u) = ~c- 1 dA(2oj)/dt. Making use of Eq. ^15j) one 
can see that the Stokes parameters of emitted radiation 
are given by 



P,(2lu) 
Pl(2u) 
P c {2lo) 



| Ja (2o;)| 2 -|jW| 2 
| J;c (2c)| 2 + | J ,(2 W )| 2 ' 
2Re{j x (2uj)f y (2Lj)} 
\j x (2u>)\* + \j y (2u>W' 
2lm{j x (2uj)f y (2uj)} 
|^(2 W )| 2 + |j,(2 W )| 2 



(16) 



Here Pi and P[ are the degrees of linear polarization in 
two coordinate frames: the xy frame and in the x'y' 
frame which is rotated by 45° with respect to xy frame. 
Quantity P c is the circular polarization degree. In deriva- 
tion of Eqs. (16) we assumed small incidence angles, 
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FIG. 3: Frequency dependence of Stokes parameters of sec- 
ond harmonic: Pi is shown by black, P{ by blue and P c by 
red curves, respectively. It is assumed that the fundamen- 
tal harmonic is fully linearly polarized in x'y' axis frame: 
\E X \ 2 = \E y \ 2 , 2Re{E x E*} = \E X \ 2 + \E y \ 2 . Calculation is 
carried out for degenerate electrons in the model of the short- 
range scattering where t\ = 2t2 oc e^ 1 [panel (a)] and of the 
Coulomb scattering where n = 3t2 oc et [panel (b)]. 



q\\ <C q z , and considered the emitted wave propagating 
along positive z axis. 

We assume that the in-plane wave vector is directed 
along x axis, q x 7^ 0, q y = 0. Non-trivial polarization 
effects take place if both j x (2u>) and j y (2uj) are non-zero. 
Such situation may arise for the linearly polarized inci- 
dent radiation in x'y' axis frame, where \E X \ 2 = |-E a | 2 , 
2Re{E x E*} = \E X \ 2 + \E y \ 2 [Pi(u) = 0, i*(w) = 1, 
P c (lu) = 0] which is considered in the following. Stokes 
parameters for the second harmonic are plotted in Fig. [3] 
as functions of the fundamental radiation frequency. Cal- 
culation presented in Fig. [3]ja) is carried out for the short- 



range scattering where t\ = 2r 2 cx e^ 1 [panel (a)]. For 
(15) comparison Fig. |3jb) shows calculation results for the 
unscreened Coulomb scattering potential, where Vfc_ p cx 



\k— in Eq. ( 10 ) and one can show that t\ = 3t 2 oc Sk- 
in all calculations it is assumed that the electron gas is 
degenerate. 

Figure [3] clearly demonstrates that all Stokes param- 
eters are non-zero, in general. In the limit of static 
fields, ujt\ = 0, cllipticity is forbidden, P c (2lu) = in 
agreement with Fig. ;3l In the limit of high frequencies 
UTx,WT2 3> 1 (but hio <C £f) constants Si and S2 are 
related as S2 — Si/ 2 which yields 



Pi 



P c = 0, at uti,ujt 2 > 1. 



It follows from Fig. [3] that second harmonic radiation 
is, in general, clliptically polarized. The remarkable de- 
gree of circular polarization appears at ojti ~ 1. For the 
short-range scattering considered in panel (a) of Fig. [3] 
the maximum absolute value of the circular polarization 
|P C | is about 90% at u>n ~ 1.78. The maximum value of 
the circular polarization degree for the Coulomb scatter- 
ing is somewhat smaller, |P C | 45% at uiTy ~ 1.31. Phys- 
ically, the circularly polarized second harmonic emission 
is a result of retardation: j x (2u>) and j y (2co) are phase 
shifted with respect to each other. This phase differ- 
ence can be significant resulting in almost fully circularly 
polarized second harmonic emission. In agreement with 
qualitative considerations and Fig.[2j the phase shift van- 
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ishes both for static and high frequency fields, where Si 
and S2 are real. Circular polarization approaches zero in 
these two limits. 

The appearance of the linear polarization in xy axis 
frame, Pi(2uj) ^ 0, and of the circular polarization, 
P c {2uj) ^ can be considered as a non-linear polarization 
conversion process. Inversion of initial linear polarization 
P({oj) — ► — 1 results in the reversal of the circular polar- 
ization, P c , and of the linear polarization in x'y' axis 
frame, P/(2w). Linear polarization in xy axis frame does 
not change its sign. 

Concluding remarks. We have developed theory of 
second harmonic generation in graphene. The analy- 
sis shows that in the simplest model of conduction and 
valence bands formed by tt orbitals of carbon atoms, 
the response at the double frequency appears with al- 
lowance for the photon wave vector transfer to electron 
system and it is described by two independent parame- 



ters. The microscopic theory is presented here for the 
classical frequency range, where the photon energy is 
much smaller than the characteristic energy of charge 
carriers in graphene. In this regime effects of other or- 
bitals are negligible and the second harmonic generation 
is fully determined by the specifics of the radiation mo- 
mentum transfer to "Dirac" electrons in graphene. 

We have analyzed the polarization of the emitted ra- 
diation. It is demonstrated, that even for linearly polar- 
ized fundamental harmonic, the second harmonic can be 
strongly circularly polarized. The degree of circular po- 
larization can be controlled by the radiation frequency, 
it can be reversed by changing the linear polarization of 
incident wave. 

Author is grateful to E.L. Ivchenko and S.A. Tarasenko 
for valuable discussions. The support by RFBR and "Dy- 
nasty" Foundation — ICFPM is gratefully acknowledged. 
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